{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# FEM Programming on Surface Mesh"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "这一节主要讨论曲面有限元程序实现中的数学细节，这些数学细节是 fealpy 中进行正确高效程序实现的基础。"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 符号"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "|Notation | Meaning |\n",
    "|:-----|:----- |\n",
    "|$S$ | $\\mathbb R^3$ 空间中的曲面 |\n",
    "|$K\\subset \\mathbb R^2$ | 二维空间中的标准单元|\n",
    "|$\\mathbf u = (u, v)^T$  | 二维空间中的坐标系 |\n",
    "|$\\tau_h \\subset R^3$ | 三维空间中的尺寸为 $h$ 的平面三角形，假设它的三个顶点在曲面 $S$ 上 |\n",
    "|$\\mathbf x = (x, y, z)^T\\in \\tau_h$ | $\\tau_h$ 上的一个点 |\n",
    "|$\\mathcal P_0$ | $S$ 邻近区域到 $S$ 的投影 |\n",
    "|$\\mathbf x_i,  i=1,\\cdots, n_{dof}$ | $\\tau_h$ 上 $p$ 次 Lagrangian 基函数对应的自由度坐标点, 假设 $x_i \\in S$\n",
    "|$\\tau_p\\subset \\mathbb R^3$ | 定义在 $\\tau_h$ 上的 $p$ 次多项式曲面三角形 |\n",
    "|$ \\mathbf x_p =(x_p, y_p, z_p)^T \\in \\tau_p$ | $\\tau_p$ 上一个点的三维坐标 |\n",
    "|$\\tau_S\\subset \\mathbb R^3$ | 把 $\\tau_h$ 投影到曲面 $S$ 上的曲面三角形 |\n",
    "|$ \\mathbf x_S =(x_S, y_S, z_S)^T \\in \\tau_S$ | $\\tau_S$ 上一个点的三维坐标 |\n",
    "|$\\varphi_i(\\mathbf x)$ | 定义在 $\\tau_h$ 上第 $i$ 个 Lagrangian 基函数 |"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "\n",
    "## $\\tau_h$, $\\tau_p$ 和 $\\tau_S$ 之间关系 "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "对于 $\\tau_p$ 上的任意一点 $\\mathbf x_p$, 存在一点 $\\mathbf x \\in \\tau_h$， 使得\n",
    "\n",
    "$$\n",
    "\\mathbf x_p = \\sum_{i=1}^{n_{dof}} \\mathbf x_i \\varphi_i(\\mathbf x)\n",
    "$$\n",
    "\n",
    "进一步， 存在标准参考单元 $K$ 中存在一点 $\\mathbf u = (u,v)$， 可得\n",
    "\n",
    "$$\n",
    "\\mathbf x(u,v) = \\lambda_0 \\mathbf x_0 + \\lambda_1 \\mathbf x_1 + \\lambda_2 \\mathbf x_2\n",
    "$$\n",
    "其中 $\\mathbf x_0$, $\\mathbf x_1$ 和 $\\mathbf x_2$ 为$\\tau_h$ 的三个顶点, \n",
    "$$\n",
    "\\lambda_0 = 1- u - v， \\lambda_1 = u, \\lambda_2 = v\n",
    "$$\n",
    "\n",
    "对于 $\\tau_S$ 上的任意一点 $\\mathbf x_S$, 存在 $\\tau_p$ 上的一点 $\\mathbf x_p$， 使得\n",
    "$$\n",
    "\\mathbf x_S = \\mathcal P_0(\\mathbf x_p)\n",
    "$$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "$\\mathbf x$ 关于  $(u, v)$ 的 Jacobi 矩阵为\n",
    "$$\n",
    "\\frac{\\partial \\mathbf x}{\\partial \\mathbf u} = [\\mathbf x_1 - \\mathbf x_0, \\mathbf x_2 - \\mathbf x_0]\n",
    "$$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "$\\mathbf x_p$ 关于 $\\mathbf x$ 的 Jacobi 矩阵为\n",
    "$$\n",
    "\\frac{\\partial \\mathbf x_p}{\\partial \\mathbf x} = \\sum_{i=1}^{n_{dof}}\n",
    "\\begin{bmatrix}\n",
    "x_i\\nabla_{\\mathbf x}\\varphi_i(\\mathbf x)^T\\\\\n",
    "y_i\\nabla_{\\mathbf x}\\varphi_i(\\mathbf x)^T\\\\\n",
    "z_i\\nabla_{\\mathbf x}\\varphi_i(\\mathbf x)^T\\\\\n",
    "\\end{bmatrix}\n",
    "$$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "则 $\\mathbf x_p$ 关于 $\\mathbf u$ 的 Jacobi 矩阵为\n",
    "$$\n",
    "\\frac{\\partial \\mathbf x_p}{\\partial \\mathbf u}=[\\frac{\\partial \\mathbf x_p}{\\partial u}, \\frac{\\partial \\mathbf x_p}{\\partial v}]=\\sum_{i=1}^{n_{dof}}\n",
    "\\begin{bmatrix}\n",
    "x_i\\nabla_{\\mathbf x}\\varphi_i(\\mathbf x)^T\\\\\n",
    "y_i\\nabla_{\\mathbf x}\\varphi_i(\\mathbf x)^T\\\\\n",
    "z_i\\nabla_{\\mathbf x}\\varphi_i(\\mathbf x)^T\\\\\n",
    "\\end{bmatrix}\n",
    "[\\mathbf x_1 - \\mathbf x_0, \\mathbf x_2 - \\mathbf x_0]\\qquad(1)\n",
    "$$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "记\n",
    "$$ \n",
    "\\mathrm d \\mathbf x_p = \\frac{\\partial \\mathbf x_p}{\\partial \\mathbf u}\\mathrm d \\mathbf u = \\frac{\\partial \\mathbf x_p}{\\partial u}\\mathrm d u + \\frac{\\partial \\mathbf x_p}{\\partial v}\\mathrm d v,\n",
    "$$\n",
    "其中 $\\mathrm d \\mathbf u = [\\mathrm d u, \\mathrm d v]^T$。\n",
    "\n",
    "进一步可得曲面三角形 $\\tau_p$ 上的第一基本形式\n",
    "$$\n",
    "I = <\\mathrm d \\mathbf x_p, \\mathrm d \\mathbf x_p> = \\mathrm d \\mathbf u^T \n",
    "\\begin{bmatrix}\n",
    "g_{11} & g_{12}\\\\\n",
    "g_{12} & g_{22}\n",
    "\\end{bmatrix}\n",
    "\\mathrm d \\mathbf u\n",
    "$$\n",
    "其中 \n",
    "$$\n",
    "g_{11} =<\\frac{\\partial \\mathbf x_p}{\\partial u}, \\frac{\\partial \\mathbf x_p}{\\partial u}>, \n",
    "g_{12} =<\\frac{\\partial \\mathbf x_p}{\\partial u}, \\frac{\\partial \\mathbf x_p}{\\partial v}>, \n",
    "g_{22} =<\\frac{\\partial \\mathbf x_p}{\\partial v}, \\frac{\\partial \\mathbf x_p}{\\partial v}>, \n",
    "$$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "定义$\\tau_p$ 上的基函数如下\n",
    "$$\n",
    "\\varphi_{p,i}(\\mathbf x_p) =\\varphi_i(\\mathbf x) \n",
    "$$\n",
    "其中\n",
    "$$\n",
    "\\mathbf x_p = \\sum_{i=1}^{n_{dof}} \\mathbf x_i \\varphi_i(\\mathbf x)\n",
    "$$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "则 $\\varphi_{p,i}(\\mathbf x_p)$ 在 $\\tau_p$ 上的切向导数定义如下：\n",
    "$$\n",
    "\\nabla_{S_p} \\varphi_{p,i} = \\frac{\\partial \\mathbf x_p}{\\partial \\mathbf u}\\begin{bmatrix}\n",
    "g_{11} & g_{12}\\\\\n",
    "g_{12} & g_{22}\n",
    "\\end{bmatrix}^{-1}(\\frac{\\partial \\mathbf x}{\\partial \\mathbf u})^T\\nabla_{S_h}\\varphi_i(\\mathbf x)\n",
    "$$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## $S$ 上曲面三角形的面积计算公式"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "$$\n",
    "\\mathcal P_0(\\mathbf x):=\\mathbf x - d(\\mathbf x)\\mathbf n(\\mathbf x)\n",
    "$$\n",
    "\n",
    "对于  $\\mathbf x_p \\in \\tau_p$, 存在 $\\mathbf x_S \\in S$, 有\n",
    "$$\n",
    "\\mathbf x_S = \\mathcal P_0(\\mathbf x_p)=\\mathbf x_p - d(\\mathbf x_p)\\mathbf n(\\mathbf x_p)\n",
    "$$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "$$\n",
    "\\frac{\\partial \\mathbf x_S}{\\partial\\mathbf x_p} = I - d(\\mathbf x_p) H(\\mathbf x_p) - \\mathbf n(\\mathbf x_p)\\mathbf n(\\mathbf x_p)^T\n",
    "$$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "$$\n",
    "\\frac{\\partial \\mathbf x_S}{\\partial\\mathbf u} = \\frac{\\partial \\mathbf x_S}{\\partial\\mathbf x_p}\\frac{\\partial \\mathbf x_p}{\\partial\\mathbf u}\n",
    "$$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## $S$ 上的导数计算"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "考虑 $\\tau_S$ 和 $\\tau_p$ 的关系\n",
    "\n",
    "则 $\\mathbf x_S$ 关于 $\\mathbf u$ 的 Jacobi 矩阵为\n",
    "\n",
    "$$\n",
    "\\frac{\\partial \\mathbf x_S}{\\partial\\mathbf u} = \\frac{\\partial \\mathbf x_S}{\\partial\\mathbf x_p}\\frac{\\partial \\mathbf x_p}{\\partial\\mathbf u}\n",
    "$$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "$$\n",
    "\\frac{\\partial \\mathbf x_S}{\\partial\\mathbf x_p} = I - d(\\mathbf x_p) H(\\mathbf x_p) - \\mathbf n(\\mathbf x_p)\\mathbf n(\\mathbf x_p)^T\n",
    "$$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "由于\n",
    "\n",
    "$$\n",
    "\\frac{\\partial \\mathbf x_p}{\\partial \\mathbf u}=[\\frac{\\partial \\mathbf x_p}{\\partial u}, \\frac{\\partial \\mathbf x_p}{\\partial v}]=\\sum_{i=1}^{n_{dof}}\n",
    "\\begin{bmatrix}\n",
    "x_i\\nabla_{\\mathbf x}\\varphi_i(\\mathbf x)^T\\\\\n",
    "y_i\\nabla_{\\mathbf x}\\varphi_i(\\mathbf x)^T\\\\\n",
    "z_i\\nabla_{\\mathbf x}\\varphi_i(\\mathbf x)^T\\\\\n",
    "\\end{bmatrix}\n",
    "[\\mathbf x_1 - \\mathbf x_0, \\mathbf x_2 - \\mathbf x_0]\n",
    "$$\n",
    "即,很容易计算\n",
    "$$\n",
    "\\frac{\\partial \\mathbf x_S}{\\partial\\mathbf u}\n",
    "$$ "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "进一步可得到曲面 $\\tau_S$ 上的第一基本形式\n",
    "\n",
    "记\n",
    "$$ \n",
    "\\mathrm d \\mathbf x_S = \\frac{\\partial \\mathbf x_S}{\\partial \\mathbf u}\\mathrm d \\mathbf u = \\frac{\\partial \\mathbf x_S}{\\partial u}\\mathrm d u + \\frac{\\partial \\mathbf x_S}{\\partial v}\\mathrm d v,\n",
    "$$\n",
    "\n",
    "其中\n",
    "\\begin{align}\n",
    "\\mathrm d \\mathbf x_S & = \\frac{\\partial \\mathbf x_S}{\\partial \\mathbf u}\\mathrm d \\mathbf u \\\\\n",
    "& =\\frac{\\partial \\mathbf x_S}{\\partial\\mathbf x_p}\\frac{\\partial \\mathbf x_p}{\\partial\\mathbf u}\\mathrm d \\mathbf u \\\\\n",
    "& = \\frac{\\partial \\mathbf x_S}{\\partial\\mathbf x_p}\\frac{\\partial \\mathbf x_p}{\\partial u}\\mathrm d u + \\frac{\\partial \\mathbf x_S}{\\partial\\mathbf x_p}\\frac{\\partial \\mathbf x_p}{\\partial v}\\mathrm d v\\\\\n",
    "\\end{align}\n",
    "故\n",
    "\\begin{align}\n",
    "\\frac{\\partial \\mathbf x_S}{\\partial u} &= \\frac{\\partial \\mathbf x_S}{\\partial\\mathbf x_p}\\frac{\\partial \\mathbf x_p}{\\partial u}\\\\\n",
    "\\frac{\\partial \\mathbf x_S}{\\partial v} &= \\frac{\\partial \\mathbf x_S}{\\partial\\mathbf x_p}\\frac{\\partial \\mathbf x_p}{\\partial v}\n",
    "\\end{align}"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "$$\n",
    "I = <\\mathrm d \\mathbf x_S, \\mathrm d \\mathbf x_S> = \\mathrm d \\mathbf u^T \n",
    "\\begin{bmatrix}\n",
    "g'_{11} & g'_{12}\\\\\n",
    "g'_{12} & g'_{22}\n",
    "\\end{bmatrix}\n",
    "\\mathrm d \\mathbf u\n",
    "$$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "其中 \n",
    "$$\n",
    "g'_{11} =<\\frac{\\partial \\mathbf x_S}{\\partial u}, \\frac{\\partial \\mathbf x_S}{\\partial u}>, \n",
    "g'_{12} =<\\frac{\\partial \\mathbf x_S}{\\partial u}, \\frac{\\partial \\mathbf x_S}{\\partial v}>, \n",
    "g'_{22} =<\\frac{\\partial \\mathbf x_S}{\\partial v}, \\frac{\\partial \\mathbf x_S}{\\partial v}>, \n",
    "$$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "定义$\\tau_S$ 上的基函数如下\n",
    "$$\n",
    "\\varphi_{S,i}(\\mathbf x_S) =\\varphi_{i}(\\mathbf x) \n",
    "$$\n",
    "其中\n",
    "$$\n",
    "\\mathbf x_S = \\sum_{i=1}^{n_{dof}} \\mathbf x_i \\varphi_{i}(\\mathbf x) \n",
    "$$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "则 $\\varphi_{S,i}(\\mathbf x_S)$ 在 $\\tau_S$ 上的导数定义如下：\n",
    "$$\n",
    "\\nabla_{S_S} \\varphi_{S,i} = \\frac{\\partial \\mathbf x_S}{\\partial \\mathbf u}\\begin{bmatrix}\n",
    "g'_{11} & g'_{12}\\\\\n",
    "g'_{12} & g'_{22}\n",
    "\\end{bmatrix}^{-1}(\\frac{\\partial \\mathbf x}{\\partial \\mathbf u})^T\\nabla_{S_h}\\varphi_{i}(\\mathbf x) \n",
    "$$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "设 $w(\\mathbf x_S)$ 是定义在 $S$ 上的函数， 利用投影可以定义 $S_p$ 上函数 \n",
    "\n",
    "$$\n",
    "\\hat w(\\mathbf x_p) = w(\\mathcal P_0(x_p))\n",
    "$$\n",
    "\n",
    "下面讨论如何计算 $\\nabla_{S_p} w$. \n",
    "$$\n",
    "\\nabla_{S_p} \\hat w(\\mathbf x_p) = \\frac{\\partial \\mathbf x_p}{\\partial \\mathbf u}\\begin{bmatrix}\n",
    "g_{11} & g_{12}\\\\\n",
    "g_{12} & g_{22}\n",
    "\\end{bmatrix}^{-1}\n",
    "\\begin{pmatrix}\n",
    "\\hat w_u \\\\ \\hat w_v\n",
    "\\end{pmatrix}\n",
    "$$\n",
    "\n",
    "$$\n",
    "\\begin{pmatrix}\n",
    "\\hat w_u \\\\ \\hat w_v\n",
    "\\end{pmatrix}\n",
    "= (\\frac{\\partial \\mathbf x_S}{\\partial \\mathbf u})^T \\nabla_{\\mathbf x_S} w(x_S)\n",
    "$$"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.6.7"
  },
  "latex_envs": {
   "bibliofile": "biblio.bib",
   "cite_by": "apalike",
   "current_citInitial": 1,
   "eqLabelWithNumbers": true,
   "eqNumInitial": 0
  }
 },
 "nbformat": 4,
 "nbformat_minor": 1
}
